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Introduction o   The  problem  studied  in  this  paper  is  that  of 
obtaining  optimum  trajectories,  such  as  the  one  for  effecting 
interception  of  a  target  following  a  known  course  above  the  sen- 
sible atmosphere,  with  minimum  fuel  consumption,  or  in  minimum 
timec   A  procedure  is  given  for  determining  the  trajectory  on  a 
digital  computer,  and  some  conditions  are  given  for  the  corres- 
ponding optimum  trajectories,.   The  methods  apply  generally  to 
problems  wherein  the  range  and  time  are  not  too  great  and  the 
terminal  velocity  is  irrelevant0 

An  essential  feature  is  in  the  use  of  the  adjoint  system  of 
differential  equations  as  defined  by  Ge  A0  Bliss „   A  somewhat 
novel  feature  is  in  the  use  of  direct  methods  and  an  optimizing 
principle  in  determining  the  trajectories,,   A  simplification  of 
the  differential  equations  is  effected  by  using  the  acceleration 
rather  than  the  mass  as  a  variable . 

The  steering  equation,  which  is  well  known,  is  a  simple  con- 
sequence of  the  optimum  principle  in  control  theory  which  has 
received  considerable  attention  in  this  country  and  in  Russia  re- 
cently,,  The  throttling  relation,  also  obtained  from  an  optimum 
principle,  is  apparently  not  generally  known:   in  problems  where 
the  energy  available  is  limited,  the  energy  input  is  to  occur 
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when  a  vector  defined  by  the  adjoint  system  of  differential 
equations  is  a  maximum  in  magnitude „   These  principles  reduce 
the  problem  to  the  solution  of  a  three-point  boundary  problem , 

The  following  simplifying  assumptions  are  made0   The  action 
takes  place  above  the  atmosphere  so  that  aerodynamic  forces  are 
negligible,,   The  trajectory  of  the  missile  is  assumedto  lie  in 
a  plane  containing  the  center  of  the  earth  so  that  the  trajectory 
is  two  dimensional,  to  reduce  the  number  of  variables „ 

1,  Basic  equations o   The  equation  of  motion  of  a  rocket  in  a 
gravitational  field,  subject  to  no  outside  forces,  may  be  writ- 
ten as 

OO 

[1]       r  =  g  +  a, 

where  r  is  the  position  vector,   g  is  the  acceleration  due  to 

gravity,   a  is  the  acceleration  due  to  thrust,  and  a  dot  (°)  over 

a  variable  indicates  its  time  derivative,,   We  may  write 

[2]       a  =  c»me/(l-m), 

where  m  is  the  ratio  of  the  mass  of  fuel  which  has  been  consumed 

to  the  initial  gross  mass  of  the  rocket      e  is  the  unit  vector 

in  the  direction  of  the  thrust^  and  c'   is  a  constant,, 

A  useful  kinematic  relation  connecting  the  acceleration  and 
the  fuel  consumption  is  the  following 

„t   o 

[3] 


'  a  dt  =  c»  f  ~-dt  x  -c'ln(l-m) 

Ho         Jo  ±~m 


Since  the  fuel  consumed  is  proportional  to  m  and  ln(l-m)   is  a 
mono tonic  function  of  m,   conditions  involving  the  final  value  of 
the  mass  may  be  rephrased  in  terms  of  the  integral  of  the  accelera* 
tion,  subject  to  the  constraints  on  the  size  of  mQ   For  practical 
purposes  we  may  take  m  as  bounded  above  by  some  constant  kl 


max 
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and  bounded  below  by  zero,  and  it  may  be  chosen  anywhere  on  or 
between  these  limits,  so  long  as  any  fuel  remains. 
The  equations  of  motion  then  have  the  form 

rx     =  gi  +  a  cos  p 

W   {.. 

y  =  g2  +  a  sin  p  , 

where  we  may  think  of  gt     as  -c  x/r3  ,  with  c   a  constant  and 

©  o 

r2  =  x2+  y2 ,  and  g2   as  -c  y/r3;  these  equations  may  be  refined 
or  simplified o 

To  get  the  needed  formulas,  let  us  multiply  equations  [4] 
through  by  two  new  variables  u,v,  which  are  unspecified  so 
far,  and  integrate  formally  to  get 

[5]   [u(x  -gx  -a  cos  p)  +  v(y  ~g2   -a  sin  p) ]dt  =  0  , 

•'o 

where  x,y  are  any  solutions  to  [4]0 

Let  us  consider  also  a  neighboring  path  whereon  the  original 

equations  are  satisfied 0   For  the  first  variations,  equation  [5] 

becomes 

T 
[6]  [u(6x  -   g!    6x  °glTr6y  -   6a  cos  p  +   a  6p  sin  p) 

J0  y 

+  v(fiy   -   g2x6x  °g2yfy  -   *a  sin  P  ■=   a  ^p  cos  p)  ]dt  =  0; 

subscripts  denote  partial  derivatives     g±      =  9gi/Sx,    etc0      If 
we  integrate  by  parts,    this  may  be  rewritten 

[7]      [u6x  -  u6x  +  v6y  -  v6yj+   [a(u  cos   p  +  v   sin  p)  ]        6t± 

T 
+     J  I>x(u  =glxu  -g2xv)    +   6yfv   -   glyu  -   g2yv) 

-  6a(u  cos  p  +  v   sin  p) 

-  a  6p(-u   sin  p  +  v   cos  p) ]dt  =  0; 

here      tx      is   a  symbol   for   any  time  or   times  when     a     is 
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discontinuous , 

Now,  to  simplify  equation  [7],  let  us  choose  u,v  as  solu- 
tions to  the  system  of  differential  equations 

r"  "  g*xu  "  g2xV  =  ° 
[8]      {00 

v  -  giyu  -  g2yv  =  0  ; 

this  choice  of  u,v  eliminates  the  variations  of  the  dependent 
variables  x,y  from  the  integral  in  [7]„   Equations  [8]  are 
called  adjoint  to  the  variational  equations  corresponding  to 
[4];  in  general,  it  is  an  integration  by  parts  and  setting  to 
zero  the  coefficients  of  terms  involving  the  dependent  varia- 
bles which  determines  the  adjoint  systemQ   Let  us  further  choose 
the  solutions  so  that  u(T)  =  v(T)  =  0o   Equation  [7]  is  then 
a  fairly  general  relation  connecting  the  variations  of  the  end 
values  of  x,y  with  the  variations  of  the  control  variations 
a,p  for  any  particular  fixed  value  T  and  any  particular  trajec- 
tory. 

2„  Variations  of  dependent  variables 0   If  we  further  choose  the 

solution  u^Vj.   to  [8]  such  that 

(u1   =  -1 
[9]      {. 

Vi  =   0, 

and  if  all  of  the  initial  values  are  given,  equation  [7]  becomes 

[10]  6x(T)  =  -[a(u±cos  p  +  ViSin  p)]£1+6t,  + 

T 
+       [SaCUiCos  p  +  vxsin  p)    +   a  £p(-Ui_sin  p  +  Vicos  p)]dt; 
J0 

similarly,   if     u2(T)    =  0,   v2(T)    =  -1, 


[11]    6y(T)    =  -[a(u2cos  p  +  v2sin  p)]:1    6tx 

t1» 


+ 


-5-  RP  29,  DS-11641 

I 

[6a(u2cos  p  +  VgSin  p)   +   a  6p(=u2sin  p  +  v2cos  p)]dt0 
t) 

These   are  the  essential  equations  for  variation  and   control, 

since  they  express   the  effect  on  the   terminal  values  of     x,y 

of   small  changes  in   the   control  variables     a,pe 

It  will  be  necessary  to  have   a  fundamental   set  of  solutions 

112      2 
to   equations   [8];   we  may  choose,    to  be  definite,    u   ,v    ,u    ,v    , 

3     3     4     4 
u    ,v    ,u   ,v  ,     such  that 

AT^O)    =  1,   v1(0)    =  0,   u^O)    =  0,   v1(0)    =  0 

u2(0)    =  0,   v2(0)    =  1,   u2(0)    -  0,   v2(0)    =  0 
[12]  <    3 

u   (0)    =  0,°°° 

u4(0)    =  0,'»  v4(0)    =  1    c 

Then,   for  any  particular   trajectory  and   any  time     T,    every 

solution  is   a  linear  combination  of  these,    and  we  may  write 

/ui  =  L  Gi^ 

[13]  J  \  i   =  1,2 

w.  =  y  c  .vj  , 

1 
for  the  proper  choice  of  the  c's0 

3,  Maximizing  principle0   It  is  a  property  of  extremals  in  nor- 
mal problems  that  they  furnish  a  maximum  to  an  integral  and  that 
this  is  done  by  maximizing  the  integrand;  this  may  be  taken  as 
the  characteristic  of  an  extremal  as  follows 0 
Let  us  write  [5]  in  the  alternate  form 

[14]   [ux  -  ux  +  vy  -  vy]n  +   (xu  -  uga  +  yv  -  vg2)dt 


'0 


0 

rT 


a°w  dt 
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Now  suppose  we  have  selected  w  =  ui  +  v"J  in  some  way,  where 

u,v   satisfy  [8]  and  w  =  |w(t)|   is  a  diminishing  function  of 

t  with  w(T)  =  Oo   Without  regard  for  a  reason,  let  us  consider 

the  problem  of  maximizing  the  integral  on  the  right  in  [14],  by 

choosing  a  properly 0   It  most  satisfy  the  constraints  that  (i) 

0  ^  a  <,   a   ,  where  a     is  a  known  function  of  t,  and  (ii) 
max7        max  * 

[15]       a  dt  ^  C,  a  given  constant; 

J0 

C:  is  by  [3]  a  measure  of  the  fuel  consumption,,   For   a     we  may 
take  the  function 

[16]      amax  =  ^"^^'W^Aax*5*  cW(1Aax)]' 

where  m  ,  *&„,_  are  the  maximum  values  of  m,m,  determined 
max'  max  '  ' 

by  the  rocket  mechanism „ 

We  see  that  to  maximize  the  integral  in  [14]  it  is  neces- 
sary and  sufficient  that  (i)   a  is  parallel  to  w,  so  that 

[17]      tan  p  =  v/u 

(the  proper  branch  of  p  chosen)  and  (ii) 
[18]      a  =  amax  ,   Octet* 
H)   ,    ta  *-  t  *  Z, 

tx     being  chosen  so  that  equation  [15]  is  satisfiedc 

An  outline  of  the  proof  is  given  here;  it  is  given  in 

detail  in  (1,  section  5) „   Consider  the  path  defined  above  and 

any  other  path;  let  the  functions  for  a,p  be  denoted  by  A,P 

on  the  second  patho   A  must  satisfy  the  constraints  (i),(ii) 

specified  on  a  above,  so  that 

T 
[19]      f  A  dt  s  Co 

It  is  seen  that 
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T  t  t 

(a-l)°w  dt  =    (a  2  -A)w  dt  +  f  A(l-cos[p-P])w  dt 


'0 


max 


0 


T 
A  w  cos(p-P)  dt  0 

The  sum  of  the  first  and  third  terms  is  positive  or  zero  and 
the  second  term  is  positive  or  zero.   This  establishes  that  the 
path  just  described  maximizes  the  integral . 

40  Attaining  a  fixed  point  in  specified  time  with  minimum  fuel. 
Let  us  consider,  as  an  example ,  the  problem  of  determining  a 
trajectory  to  attain  a  specified  point  X~,Y~  in  a  specified 
time  T  with  a  minimum  of  fuel,  assuming  that  all  of  the  initial 
conditions  are  givenQ   We  will  first  take  up  a  routine  for 
finding  the  trajectory  and  then  show  how  it  satisfies  certain 
conditions  which  can  be  checked 0 

Let  us  guess  initial  values  for  the  solutions  u^Vi   and 
u2»v2>  solutions  to  [8],  and  a  linear  combination  of  these 
,u  =  UiCos  ©  +  u2  sin  © 


[20] 

v  =  ViCos  9  +  v2  sin  ©, 

where  ©  is  a  number  to  be  found „   Each  of  the  quantities  should 
have  an  iteration  number,  since  they  will  be  in  an  iterative  com- 
putational routine „  Now  let 
[17]      tan  p  =  u/v, 

properly  chosen,,   Let  us  guess  also  t^  maximum  thrust  is  to  be 
applied  until  time  tl9    and  thrust  is  zero  thereafter „ 

Now  compute  the  corresponding  trajectory,.   It  will  lead  to 
terminal  values  x(T),y(T),  which  are  in  error  by  amounts 


Xf-x(T) ,Yf»y(T) c   Compute  simultaneously  a  fundamental  set  of 
solutions  u  ,v   to  the  adjoint  system  of  differential  equations 
[8],  and  two  further  integrals  which  occur  later  in  equation  [22] 
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From  [17], [20],  for  fixed  or  determined  u1,v1,u2,v2,  we  get 

[21]      6p  =  bQiu^z   =  u2v1)/([u]2+[v]2)0 

The  relations  for  the  variations  of  the  end  values  then  become 


6x(T)  =  [a(u!Cos  p  +  v±sin  p) ]^  _6tj. 

[(-^sin  p  +  Vicos  p)(uiv,-u9vi)dt  6G 
[u]2+  [v]2 


ai 
0 


[22],  . 

^6y(T)  =  [a(u2cos  p  +  v2sin  p)]^  _6ti  +  J, 

Now  set 

r6x(T)  =  X,  -  x(T) 

[23]     {         f 

L6y(T)  =  Yf  -  y(T)  , 

substitute  into  [22],  and  solve  for  6tl76©.   This  gives  corrected 

values  of  t1?9  for  the  next  round0 

1     4 
From  the  fundamental  set  of  solutions  u  ,  °°°v  ,   corrected 

estimates  are  made  for  the  initial  values  of  'Ui?v1,u2,v2   and 

their  derivatives „   All  of  the  initial  values  are  then  available 

for  starting  the  next  iterationc 

This  routine  is  continued  until  some  convergence  criterion 
is  satisfied:   for  example  one  might  require 
[24]      E  =  (Xf  -  x[T])2  +(Yf  -  y[T])2 
to  be  less  than  some  preassigned  number „ 

If  the  routine  given  above  converges,  we  have  found  a 
trajectory  which  effects  at  least  a  stationary  value  for  the 
fuel  consumption  compared  with  paths  in  some  neighborhood 0   For 
if  we  consider  neighboring  paths  in  equation  [14],  we  see  that 
the  first  variation  of  the  left  side  vanishes,  by  virtue  of 
equations  [8]  and  the  fact  that  6x(0) ,6y(0) ,6x(0) ,6y(0) ,6x(T) , 
6y(T) ,u(T) ,v(T)   all  vanisho  But  the  right  side  of  the  equation 
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is  positive  for  every  admissible  variation  of  the  acceleration, 
as  was  shown  in  section  3C   Hence  there  is  no  neighboring  path 
whose  first  variations  do  not  vanish  which  requires  the  same 
or  less  fuel.   If  gravity  can  be  approximated  as  a  linear  func- 
tion of  position,  this  also  establishes  the  sufficient  conditions 
for  an  absolute  minimum  [see  (1),  section  5]0 

5.  Variations  of  terminal  time;  minimum  fuel  consumption 0  Now 
suppose  that  the  terminal  time  T  may  vary0   The  change  in  the 
end  values  of  x,y  are 

r6X  =  x(T)6T  +  6x(T) 
[25]     { 

L6Y  =  y(T)6T  +  6y(T)  . 

There  are  corresponding  equations  for  the  terminal  values  of 

u1?  etCo, 

Ui(T)  =  -UiCT^T,  Vl(T)  =  -Vl(T)6T,  u±(T)  =  -1  -ui(T)6T, 

etc;  in  these  equations  T  should  have  an  iteration  subscript, 

being  the  approximation  associated  with  the  last  trajectory  run, 

Since  there  are  three  variables  t1?T,9,  and  only  two  equa- 
tions for  Xf,Yf,  a  further  relation  is  required,,  If  we  are  in- 
terested in  minimum  fuel  consumption,  then  tt  must  be  a  minimum . 
Equation  [7]  may  be  rewritten 

[26]      (=u6x  -  v6y)T  =  [a(u  cos  p  +  v  sin  p)].   6t±  , 
if  terms  that  must  vanish  are  omittedo 

If  we  consider  two  neighboring  trajectories,  each  of  which 
effects  interception  of  a  target  whose  coordinates  are  Xf(T), 
Yf(T),  then  6x(T)  =  (Xf~x)6T,  6y(T)  =  (Yf-y)6T,  and  [26] 
becomes 
[27]      (-w)T°[(Xrx)?  4-  (Yf-y)J]T6T  =  (aw)t  _6tlQ 
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Hence,  for  minimum  fuel  consumption, 

[28]      w(T)°[(Xf-x)?  +  (Yf-;yO^T  =  0: 

these  two  vectors  must  be  perpendicular;  apparently  this  is 
a  fairly  general  condition  associated  with  minimum  fuel  con- 
sumption,,  This  condition  allows  a  direct  determination  of  the 
angle  0O   This  in  turn  defines  69  as  the  difference  of  the 
two  successive  values  of  0O 

The  value  for  6x(T)  which  must  be  used  in  the  iterative 
routine  is  6x(T)  =  Xf(T)  -  x(T)  +  (Xf  -  x)T6T,  and  the  equa- 
tions for  corrections  to  T,tA   are 

Xf(T)  -  x(T)  =  -  (Xf  -  x)6T  +  6x(T) 


[29] 

-'Yf(T)  -  y(T)  =  -  (Yf  -  y)6T  +  6y(T)  , 

where  6x(T),6y(T)?  linear  in  6G96tl9   are  given  by  equations 
[22] o   Since   69  is  already  determined,  these  equations  [29] 
then  yield  6T,6tl0 

The  iteration  is  continued  until  some  convergence  criterion 
is  satisfied0   It  must  also  be  checked  that  w  is  a  decreasing 
function  of  time;  actually  all  that  is  required  is  that 
w(t)  >  w(tx)  when  t  <  tx   and  w(t)  <  \f(tx)    when  t  >  tl0   It 
should  be  pointed  out  that  the  relation  involving  the  end  con- 
dition, equation  [28], is  only  a  necessary  condition,.   There  do 
not  seem  to  be  any  simple  ways  to  ensure  that  a  minimum,  rather 
than  a  maximum  is  obtained;  in  the  case  of  constant  gravity  and 
a  ballistic  target,  the  trajectory  always  furnishes  a  minimum,, 

6.  Comments 0   The  procedure  seems  to  furnish  a  rather  simple  ex- 
tension of  the  procedures  which  Bliss  (3)  introduced  in  Ballistics, 
for  calculating  differentials 0   There  seems  to  be  only  one 
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peculiarity  of  the  trajectories,  namely  that  the  vector  w 

is  decreasing  in  magnitude 0 

If  a  minimum  time  of  interception  is  desired,  then  set 

tt  =  lesser(t«    ,T)0  Equations  [29]  then  furnish  the  changes 
1  -"•max*      ^ 

6T,60o   It  is  necessary  to  make  the  first  estimate  of  0  in  a 
certain  range,  else  the  routine  may  not  converge,  or        may 
converge  to  the  trajectory  which  takes  maximum  time  (l,p  17) 0 

The  use  of  the  acceleration  rather  than  the  mass  as  a  varia- 
ble simplifies  the  differential  equations ,  since  the  mass  enters 
non  linearly  in  the  differential  equations „   On  the  other  hand, 
the  bound  on  a  is  a  function  of  time,  not  a  constant,  in  the 
above  example,  and  in  more  general  cases  it  does  not  seem  possible 
to  express  the  bounds  explicitly  as  functions  of  timeQ 

There  seem  to  be  two  general  principles  associated  with 
optimum  trajectories  where  there  is  a  limited  amount  of  energy 
available,  so  that  throttling  must  occur  in  some  form  or  another 0 
The  first  of  these  is  that  the  orientation  of  the  energy-input 
vector,  the  acceleration  in  this  case,  is  such  as  to  maximize 
its  projection  onto  the  adjoint  vectorQ   The  second  is  that  the 
energy  is  to  be  expended  at  a  time  when  the  adjoint  vector  is 
large o   The  first  condition  above  leads  to  the  well-known  steering 
criterion,  apparently  first  published  by  Lawden  [see  (2)  for  refer- 
ences] 0   The  second  is  apparently  equivalent  to  the  Weierstrass 
condition  in  some  problems,  though  as  expressed  above,  it  has  no 
analog  in  classic  calculus  of  variationsc 

The  method  of  solution,  using  procedures  due  to  Bliss  (3), (4), 
(5),  has  the  nice  feature  that,  should  the  actual  trajectory 
deviate  in  some  way  from  the  planned  trajectory,  the  same  equa- 
tions for  the  variations  can  be  used  to  determine  a  corrective 
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thrust  schedule,.   The  vectors  Wj.,^  may  be  considered  the 
impulse  response  vectors  for  x(T) ,y(T)   respectively,  since 
the  dot  product  e°w.   gives  the  change  in  the  terminal  value 
of  the  corresponding  variable  due  to  a  unit  impulse  in  the 
direction  of  e  0 

In  the  computing  which  has  been  done,  a  modified  Runge- 
Kutta  routine  was  used  which  has  the  desirable  feature  that  it 
is  easy  to  change  the  time  intervals,  which  must  be  done  in  the 
neighborhood  of  tl9T„  No  convergence  problems  have  been  en~ 
countered  in  the  ballistic  missile  interception  problem0   On 
the  other  hand  it  is  apparently  necessary  to  make  a  close  initial 
estimate  in  a  corresponding  orbital  transfer  problem0   The  values 
obtained  from  the  case  where  g  varies  at  most  linearly  should 
furnish  reasonble  starting  values5  should  convergence  be  a  problem0 

Considerable  attention  has  been  attracted  recently  by  the 
maximum  principle,  due  largely  to  papers  by  Pontriagin  (6)0   To 
the  author  this  seems  to  be  the  fundamental  way  to  approach  prob- 
lems in  the  calculus  of  variations;       the  Euler  equations 
are  derived       by  effecting  a  stationary  value  for  an  integral 
through  the  e&  processes  of  calculus0   In  most  problems  there 
seems  to  be  no  practical  new  information  in  the  maximum  prin- 
ciple,,  Superficially  the  system  of  differential  equations  is 
of  lower  order  when  the  maximum  principle  is  taken  as  the  basic 
concept 0   In  the  numerical  solution 9  the  last  Euler  equation  [17] 
may  be  solved  by  Newton's  method,  which  depends  on  derivatives 
so  that  it  is  equivalent  to  an  additional  differential  equation . 

Several  authors  have  dealt-  with  the  mass  and  handled  the 
bounds  on  m  by  introducing  a  new  variable,  such  as  (|)  in  the 
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equation 

m  =  mTriov  sin2(t)o 
max     T 

In  simple  problems ,  both  the  mathematical  theory  and  the  pro- 
gramming of  the  numerical  solution  are  simpler  in  terms  of  the 
original  variables;  the  logical  decisions  are  of  a  type  built 
into  all  general-purpose  digital  computers 0 

It  seems  to  the  author  that  one  of  the  most  important 
features  in  recent  control  theory  rests  on  the  use  and  inter- 
pretation of  the  adjoint  system  of  differential  equations 0 
As  remarked  earlier,  the  foundations  were  laid  by  Bliss,  with 
more  recent  applications  and  contributions  by  Drenick  (7) , 
Tsien  (8),  Breakwell  (.9),  and  Tyndall  (10),  among  others0   The 
mathematical  contribution  to  calculus  of  variations  is  in  the 
numerical  procedures  for  solving  two- (and  more)  point  boundary- 
value  problems  on  high-speed  computers ,   Another  method  of  solu- 
tion, the  gradient  method  of  Kelley  (11)  also  utilizes  the 
adjoint  systemc   It  is  sometimes  stated  that  the  adjoint  equa- 
tions [8]  are  necessary  conditions  for  an  extremalo   This  seems 
a  logical  error:   the  variables  u,v  are  introduced  by  the 
mathematician,  not  implicit  in  the  problem,,  and  in  turn  are 
chosen  to  satisfy  [8],  because  that  choice  simplifies  the  inte- 
gral in  equation  [7]Q  Drc  SD  Ross  usually  brings  out  this  facto 
In  problems  such  as  the  above 9    the  construction  of  the  functional 
of  Pontriagin  [(6),  p  16]  is  equivalent  to  solving  the  three- 
point  problem  involving  the  original  system  of  equations  and 
the  adjoint  system '„ 
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